Full characterization of the quantum linear-zigzag transition in atomic chains 
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A string of repulsively interacting particles exhibits a phase transition to a zigzag structure, by 
reducing the transverse trap potential or the interparticle distance. The transition is driven by 
transverse, short wavelength vibrational modes. Based on the emergent symmetry Z2 it has been 
argued that this instability is a quantum phase transition, which can be mapped to an Ising model 
in transverse field. We perform an extensive Density Matrix Renormalization Group analysis of 
the behaviour at criticality and evaluate the critical exponents and the central charge with high 
precision. We thus provide strong numerical evidence confirming that the quantum linear-zigzag 
transition belongs to the critical Ising model universality class. These results show that structural 
instabilities of one-dimensional interacting atomic arrays can simulate quantum critical phenomena 
typical of ferromagnetic systems. 
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Quantum simulators [I] are acquiring increasing 
prominence in the area of quantum technologies. The 
basic idea, dating back to Feynman [2 , is to create a 
controllable system whose dynamics reproduces that of 
a many-body quantum mechanical model that is both 
important and very difficult to solve. Experimental real- 
isation in systems with high level of control, such as ultra 
cold atoms and ions, allows the study of the properties 
of the simulated model in a way otherwise very challeng- 
ing for classical simulations [3HS]- One prominent ex- 
ample is the simulation of quantum critical phenomena 
described by (j> 4 type of models [7HTU] , which are encoun- 
tered in solid state and high-energy models and whose 
predictions need verification [TlJ [T2] . In this context, the 
zigzag instability of interacting atomic chains is a real- 
ization of </> 4 model which can be realized in laboratory 
[T5HTB] . This transition is sketched in Fig. [I] and has 
been observed in systems of singly-charged ions confined 
by external potentials, where the instability is controlled 
either by lowering the transverse potential or the linear 
density [T71IIE]- Due to its universal properties, it can be 
observed in arrays of other repulsively interacting (quasi) 
particles [TH1 [20] , such as electrons in quantum wires [5T] , 
ultracold gases in optical lattices mutually repelling via 
the dipolar interaction [22] E3] or via the off-resonant 
coupling with a Rydberg excitation 24 , or vortices in 
quantum gases [25 . In these systems ultralow tempera- 
tures are typically achieved, so that quantum effects at 
criticality are dominant. It has been argued that the 
linear-zigzag structural instability is a quantum phase 
transition |19j . which in two dimensions can be mapped 
to an Ising model in the transverse field, describing a 
ferromagnetic transition at zero temperature [12] . This 
mapping was first proposed for Wigner crystals of elec- 
trons in quantum wires [21j . and then put forward in 




Figure 1: (color online) Linear-zigzag instability in a chain 
of interacting atoms: (left) linear and (right) zigzag config- 
uration. The transition is either controlled by changing the 
longitudinal lattice spacing a or the frequency lot of the trans- 
verse harmonic confinement |13| . 

Ref. [19] using conformal-field-theory considerations mo- 
tivated by the emergent Z 2 symmetry. These studies call 
for a precise numerical verification ruling out other pos- 
sible effects, which cannot be systematically accounted 
for when performing the mapping with analytic tools. 

In this Letter we perform an extensive Density Ma- 
trix Renormalization Group (DMRG) study [2"BH2"8] and 
demonstrate that the linear-zigzag instability belongs to 
the universality class of the Ising model in a transverse 
field. In particular, we quantify the quantum correc- 
tions to the classical linear-zigzag transition and relate 
their magnitude to experimental parameters. We remark 
that the ID lattice <j> model has been previously numer- 
ically investigated by means of other numerical meth- 
ods [HI ED], including DMRG 0T]. Our work provides 
high precision values of the critical exponents and of the 
central charge of the system, giving for the first time ir- 
refutable evidence that the discrete, quantum <f> A model 
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in one dimension belongs to the universality class of the 
Ising model with transverse field. The computed critical 
exponents are summarized in the following table together 
with the values predicted |12j : 



is a scalar with the dimension of an energy, while H reads 



Quantity 


Computed 


Theory fjJF 




Anomalous dimension 


0.258 ±0.012 


0.25 


(3 


Spont. magnetization 


0.126 ±0.011 


0.125 


V 


Correlation length 


1.03 ±0.05 


1 


c 


Central charge 


0.487 ±0.015 


0.5 



Model - DMRG program simulates the quantum criti- 
cal behaviour of a chain of repulsively interacting atoms 
with mass M, in presence of an anisotropic potential con- 
fining the motion in the x — y plane, by diagonalizing the 
Hamiltonian describing the dynamics. The atoms are 
assumed to form a regular array along the x axis com- 
posed by L sites with intcrparticlc distance a. Quan- 
tum degeneracy effects are negligible due to either the 
strength of the interaction (for ions) or the external con- 
fining potential (e.g., a deep optical lattice for ultracold 
neutral atoms) [32 . The atoms position and canonically 
conjugated momentum are (xj,yj) and (Px,j,Py,j), with 
j = 1, . . . ,L, such that [Xj,p x j] = [y 3 ,p y ,i] = ihSjj. The 
Hamiltonian reads 



Pl,i+Pl,i Mu4 2 



2M 



(1) 



where lot is the frequency of the harmonic potential in 
the y-direction, Ve(x) denotes the longitudinal confine- 
ment with a quasi-flat bulk and hard edges, C- m t gives 
the strength of the repulsive interaction. For trapped 
ions with same charge Q, a = 1, C mt — Q 2 /(47re ), and 
a is of the order of micrometers [T71 [TS] . For ultracold 
dipolar gases with dipolar moment aligned perpendicu- 
lar to the plane, a = 3 and C- ln t = Cdd/{4ir) with Cdd 
the dipolar interaction strength [23], and the interparti- 
cle distance is determined by the periodicity of an optical 
lattice along x that traps the particles deep in the Mott- 
insulator phase [3J [33] ■ 

The diagonalization of Hamiltonian is a problem of 
high numerical complexity. Since the analysis is focused 
on the quantum ground-state properties and on the low 
energy excitations close to the zigzag instability, one can 
map Eq. ([!]) to an effective one-dimensional Hamiltonian 
H e g, describing a short-range theory for the transverse 
motion [I3J [Tj|l [201 EH] • It is convenient to report H e g in 
the dimensionless form H = H e g/£o, where £q = Ci n t/a a 
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and is a discretized, anti-ferromagnetic version of the 
well-known 4 -field-theory model [30] [3T] . Now, all phys- 
ical quantities are dimensionless: y = y/a, p = p/y/MEo, 
and u = ujt / \J £o / {M a 2 ) . The values Mj^i^.z are con- 
stants, their generic dependence as a function of a and 
Ci n t is given in Ref. [20]. For instance, for ions M,\ = 
I C(3), M 2 = In 2, and M 3 = f C(5), while for dipoles 
M[ = f C(5), M> 2 - | C(3) and M> 3 = C(7), where 
((n) is the Euler-Riemann zeta function. The rescaled 
space- momentum commutator reads [yj , pi] = igSj t £ with 



n 2 



fj 



Ma 2 £ 



(3) 



The parameter g replaces the Planck constant h in the 
definition of Hamiltonian Hamiltonian H is ob- 

tained from Eq. in the thermodynamical limit, keep- 
ing a constant, and assuming that the transverse dis- 
placements are smaller than the typical lattice constant, 
y/ (y 2 ) <§C a. In this regime, close to the zigzag insta- 
bility the transverse motion is effectively decoupled from 
the longitudinal excitations [T31 [33]. Moreover, in this 
limit the low frequency part of the dispersion relation of 
Hamiltonian Eq. ([l]) is equivalent to the one ofEq. §, 
as shown in Refs. [T3j US ED]. 

The rescaled Hamiltonian description of Eq. ([2| shows 
that the emergent short-range theory depends only on 
two parameters: the rescaled transverse trap frequency Co 
and the the "effective Planck constant" g. The latter pa- 
rameter is the square root of the ratio of the kinetic over 
the interaction energy, and thus it measures the strength 
of quantum fluctuations at criticality. Landau theory is 
found when g = 0: in this limit the critical value of 
the transverse frequency lo c (0) is given by the equation 
^c(O) 2 = Mi, for which the quadratic on-site potential 
in Eq. ^ vanishes [T3J [3S] . When quantum fluctuations 
are relevant, tunneling between the two wells shifts the 
critical frequency to smaller values Co c (g) < w c (0). This 
shift has been estimated in Ref. [HI [2U] . In this work, 
among other results, we provide a precise numerical de- 
termination of L~b c . 

To perform the DMRG simulations we introduce a lo- 
cal (finite) basis at every lattice site to reduce the contin- 
uous local variable yj into a discrete algebra suitable to fit 
into DMRG analysis |28j . We solve the local part of the 
Hamiltonian H in Eq. ([2|, which is a size-independent, 
homogeneous problem. We numerically diagonalize the 
local Hamiltonian, defining an (infinite) set of eigenfunc- 
tions and we truncate the local bases retaining only the 
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Figure 2: (color online) Phase diagram of an array of ions 
in the (g,Cj 2 ) parameter space; the red (green) region locates 
the linear (zigzag) phase. Inset: displacement of the critical 
square frequency |Aa}^| as a function of g. The magenta line 
is the power-law fit described by Eq. (JSJ): |A£?c| = Ag r ', with 
r = 0.82 ± 0.001 and A = 21.9 ± 0.1. The phase diagram for 
dipoles, for instance, can be obtained v ia the rescaling Cj' 2 = 
M[ + {Cj 2 - Mi)M' 2 /M 2 and g' = gy/Mj M' 2 3 /(M 3 2 M'£), 
see text. 



Figure 3: (color online) Order parameter £ as a function of Cj, 
for different L — 100, 120 . . . 300. Bottom-left inset: rescaled 
data according to £.z,(w), characterizing f(x) (see text). Top- 
right inset: two points correlation {yjijj+e) (red dots) as a 
function of the distance £ (the black line shows a quasi-critical 
fit, see text). 



square root of the structure factor density: 



d lowest-energy states, which are used to expand the full 
many-body Hamiltonian. We check a posteriori the valid- 
ity of the approximation we introduce by keeping track of 
the populations of the single-site reduced density matri- 
ces. Typically, the matrix diagonal elements decay expo- 
nentially fast with the level index, dropping below 10~ 5 
already at d ~ 8. We therefore employ a local basis di- 
mension around 10 < d < 20. We also typically adopt 
a DMRG bondlink dimension up to m ~ 30 so that er- 
rors are kept well under control. In this setup, we handle 
simulations up to L ~ 3000 sites with open boundary 
conditions. 

Results - We first characterize the phase diagram of 
the linear-zigzag transition at the thermodynamic limit 
in the Cj and g plane. The textbook order parameter for 
distinguishing between these phases is the antiferromag- 
netic order parameter £ = lim£_ i . 00 L _1 53j=i( — ^Yvj- 
We have only direct access to finite-size systems where we 
can measure the single site displacement (jjj) and then 
perform the limit. However, as there can be no sponta- 
neous symmetry breaking in finite-size samples, (jjj) = 
for every system size L. Indeed, the two ordered-phase 
configurations interfere and the system ground state is 
their even superposition. To overcome this problem we 
compute the order parameter indirectly from the two- 
point correlations (yjije), which are insensitive to the 
symmetry breaking. More precisely, we compute the 
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where l^gQ 



is the ground state at (g,Q) and size L. 
Due to the sub-extensivity of entanglement-based quan- 
tum correlations, this definition coincides with the pre- 
vious one at the thermodynamical limit. The phase di- 
agram is determined by extrapolating the order param- 
eter in Eq. Q at the thermodynamical limit: we esti- 
mate £oo(<?)£0 = lmiL-^oo £l(<7, &) in order to discrimi- 
nate whether the (g, Cj) point lies in the linear (^ = 0) 
or in the zigzag phase > 0). Figure [2] shows the re- 
sulting phase diagram. The critical boundary appears to 
follow a power-law behavior, a numerical fit of the dis- 
placement from the classical critical point results in the 
formula 



ACj 2 c =Cj 2 c (g) 



(0) = -(21.9 ± 0.1) -g' 



0.82±0.01 



(5) 



To further characterize the quantum phase transition 
we numerically evaluate the critical exponents. In par- 
ticular, the transverse displacement y plays the role of 
the direction of the spontaneous magnetization in the 
Ising model, thus the decay rates of correlators (yjijj+i) 
as a function of the distance £ reveal the anomalous di- 
mension critical exponent 77. To extract this quantity, 
we fit (yjyj+i) far from the boundaries using equation 
{yjVj+i) — G(i) = a£~ v exp(—£/X). The exponential 
correction to the power-law decay compensates for not 
being exactly at the critical point by introducing a fi- 
nite correlation length A [55]. As it can be seen in the 
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Figure 4: Inset: Von Neumann entropy 5vn(p<) (black dots) 
as a function of the partition size I fitted via c, c' using 
equation |6| (green line). Main: c values as a function of 
L = 100, 120 .. . 300 for different cj = 0.383, 0.384 . . . 0.390 
(top to bottom) at g = 0.12 for various total lengths L and 
trap frequencies to (here at g — 0.12). At the critical point 
Cjc (red data set), c is constant with L and it estimates the 
central charge. 

topmost inset of Fig. |3j the agreement is almost perfect 
over lengths of hundreds of lattice sites. We averaged 
the fitted exponent rj for different values of g while being 
as close as possible to the phase boundary, i.e. maximiz- 
ing A. The final result is 7/ av ~ 0.258 ± 0.012, in good 
agreement with the predicted value of 1/4 [TJ] . 

The spontaneous magnetization /? and correlation 
length divergence v are extracted through a finite- 
size scaling of the order parameter defined in 
Eq. Q. Following a procedure based on renormal- 
ization group analysis, we find that near the transi- 
tion the order parameter follows the scaling ~ 
j^-P/v j _ . £i/") ; for some non-universal func- 
tion f(x) which may be model dependent [37]. It is then 
possible to plot L 7i £l as a function of (Q — w c ) L 12 , 
and then tune the exponents 71 and 72 to collapse all 
the curves onto one another. A typical example of 
such analysis is depicted in Fig. || As a result of this 
rescaling we obtain /3 cxt = 71/72 — 0.126 ± 0.011 and 
Z4xt = I/72 — 1-03 ± 0.05 in perfect agreement with the 
theory [12], which predicts 1/8 and 1 respectively. This 
finite-size scaling procedure gives also another estimate 
of the critical boundary, in perfect agreement with the 
results showed in Fig. [2] 

Finally, we measure the central charge c of the 
model 38 . It is actually possible to extract the central 
charge of a critical ground state by computing the entan- 
glement related to a left-right {1.1 «-» t + 1..L) system 
partition |39j . That is, we evaluate the Von Neumann 
Entropy <Svn of the reduced density matrix pi of either 
partition: Sy-^(pi) = — Tr [pi hipi]. It is straightforward 



to determine the partition entropy when the quantum 
state is obtained via DMRG [27] li6]. In fact, DMRG al- 
gorithms provide directly the Schmidt coefficients fJ. p (£) 
for every partition £ 6 {I..L}, from which the Von Neu- 
mann entropy Svn(w) = -£pMpC0 ln MpW is easily 
computed. The profile Sy^(pi) at the critical point grows 
as 

Svn(pi) = g log I L ■ sin — J + c , (6) 

which allows to directly fit the central charge of the sys- 
tem [3TJ U2] ■ Figure [4] displays a typical result of this 
analysis: for every system size L and trapping frequency 
Co we fit the central charge c by using Eq. Q (inset). 
We then plot the fitted values of the central charge as a 
function of the system size L and extract the value that 
is size- independent, and thus critical. We average over 
different values of g and obtain c av — 0.487±0.015, which 
is in good agreement with the predicted value of 1/2 [T2] , 

In experiments the quantum critical behaviour can be 
measured, for instance, via the structure form factor [20] 
or by quenches across the critical point [33] • The quan- 
tum disordered region is accessed for temperatures T that 
are smaller than the energy scale of the gap of the Ising 
model. For detailed discussions we refer the reader to 
Ref. [19] . Our calculation allows us to precisely deter- 
mine the frequency width of the disordered phase. This 
reads Au c = y£o / {Ma 2 ) Auo c , where ACo c is given in 
Eq. |5| as a function of the parameter g. For ions g 
is proportional to the linear density 1/a, which is lim- 
ited by the Coulomb repulsion, and takes values between 
10~ 5 and 10~ 4 giving a small critical region [19] . For 
dipolar gases and Rydberg-dressed gases the situation 
can be quite different [53] [UJ. For dipoles, for instance, 
.9 = \/ a / r o, where ?- = MCi nt /h 2 is the characteristic 
length of quantum coherence, ro < lOO^m, while long- 
range order can be realized by means of an optical lattice 
(with the molecules deep in the Mott-insulator phase, so 
that quantum fluctuations along the array can be ne- 
glected [33J) so that g > 0.01. Note that, in the absence 
of an external periodic lattice, the linear-zigzag transi- 
tion in quantum gases can be observed at rather large 
densities [231134], 

To resume, we implemented a DMRG program which 
allowed us to characterize the ground state of a one- 
dimensional array of interacting atoms close to the linear- 
zigzag transition. The universality class of the criticality 
was identified by extrapolating the critical exponents as 
well as the central charge of the quantum phase tran- 
sition. The excellent match between the predicted val- 
ues and the results of the simulations demonstrates the 
correspondence between linear-zigzag instability and the 
Ising universality class. These results show that one- 
dimensional quantum critical phenomena typical of fer- 
romagnetic systems, recently simulated by means of en- 
gineered coupling between internal and external degrees 
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of freedom of ions [5J |S] , can be naturally simulated by 
structural instabilities of interacting atomic arrays. 
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